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Abstract 

Through a suitable expansion of the Gross-Pitaevskii equation near the classi- 
cal turning point, we obtain an explicit solution for the order parameter at the 
boundary of a trapped Bose gas interacting with repulsive forces. The kinetic 
energy of the system, in terms of the classical radius R and of the harmonic 
oscillator length a HO , follows the law E^i n jN oc R~ 2 [log(R/a HO ) + const.], 
approaching, for large R, the results obtained by solving numerically the 
Gross-Pitaevskii equation. The occurrence of a Josephson-type current in the 

presence of a double trap potential is finally discussed. 
PACS numbers: 03.75.Fi, 05.30.Jp, 32.80.Pj 
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The recent experimental realization of Bose-Einstein condensation in atomic gases con- 
fined in magnetic traps [|1|-|3|] is stimulating a novel interest in the study of inhomogeneous 
Bose condensed systems where the order parameter exhibits an important spatial depen- 
dence on a macroscopic scale [[|. 

The purpose of the present work is to investigate the behavior of the wave function of 
the condensate near the classical turning point, that is, at the boundary of the trapped gas. 
This region is particularly important for the determination of the kinetic energy associated 
with the atoms of the condensate |5|,§]. It is also crucial for the description of Josephson-type 
effects taking place in the presence of a barrier in the confining potential. 

The order parameter ^(r) associated with the ground state of a dilute Bose gas obeys 
the Gross-Pitaevskii equation: 



-^V 2 + K, 4 (r) + ^hKr)| 2 
2m m 



^(r) = /#(r) , (1) 



where V ex t is the external confining potential, fi is the chemical potential and a is the s- 
wave scattering length. The condensate wave function ip(r) is normalized to the number 
iV of atoms and is related to the atomic density through p(r) = \ip(r)\ 2 . The solution of 
(|I|) has been recently found by direct numerical integration in the case of both isotropic 
and anisotropic traps ||. In the following we will consider systems interacting with 
repulsive forces (a > 0). When the scattering length (or the number of atoms in the trap) is 
sufficiently large, the solution of Eq. ([]]), in the region where \i > V ext (r), takes the simplified 
Thomas- Fermi form 



^(r) 
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Equation (fj) is obtained by neglecting the kinetic energy term V 2 ip(r) in the Schrodinger- 
like equation (|T]) and provides an accurate description of the exact solution in the interior 
of the atomic cloud where the gradients of the wave function are small, as shown in Fig. 1. 

Near the boundary region, where V ex t(r) ~ /i, the kinetic energy term in Eq. Q can not 
be longer ignored and the Thomas-Fermi approximation (0) is inadequate. This can have 



important consequences in the determination of relevant physical quantities. For example if 
one evaluates the kinetic energy associated with the condensate 

r h 2 

E kin = Jdr — |W>(r)| 2 (3) 

using Eq. (0), one finds a logarithmic divergency []5| which reveals that the boundary region 
must be investigated with higher accuracy. In the following we will explore the correct 
behavior of the order parameter in this region starting from the Gross- Pit aevskii equation 
(11). With respect to similar procedures used in the study of the single-particle Schrodinger 
equation in the presence of an external field U, the present method includes explicitly the 
interatomic forces which are responsible for crucial non linear effects in the equations of 
motion. 

Let us consider for simplicity a spherical trap || and let R be the boundary of the system, 
determined by the equation /i = V ex t(R). Near this point one can carry out the expansion 

V ext {r) -n=(r-R)F + o(r - R) (4) 

where F is the modulus of the attractive external force F = —W ext evaluated at r = R. 
Close to the boundary, where \r — R\ the Gross-Pitaevskii equation takes the form 

h 2 d 2 , , 47rfi 2 a , , , 

Let us now introduce the dimensionless variable 

where 

is a typical thickness of the boundary giving, as we will see later, the distance from the 
classical radius R where the Thomas-Fermi approximation starts failing. Then we introduce 
the adimensional function defined by 




(8) 



In terms of the Gross-Pitaevskii equation (|5|) takes the universal form 



0"-(£ + 2 )0 = O. 



(9) 



Notice that the non linear term 3 arises from the internal potential energy in the Gross- 
Pitaevskii Eq. (|l|). When £ — > +oo this term can be neglected and Eq. (||) takes the simpler 
form 0" — £0 = which is the equation defining the Airy function. The asymptotic behavior 
has then the form 



where the constant A must be determined by numerical integration of Eq. (|9]). In the opposite 
limit £ — > — oo one can neglect the second derivative 0" and the asymptotic behavior is given 



The full behavior of the function is shown in Fig. 2. The value of the constant A is found 
to be 0.397. 

The solution of Eq. (|9]) provides, via Eqs. the proper structure of the wave function 

of the condensate near the classical turning point R. It is worth noting that Eq. (§) does 
not depend on the form of the external potential nor on the size of the interatomic force. 
These physical parameters enter the transformations (|6]) and @ which fix, together with 
the solution of @, the actual behavior of the wave function ip. 

Equations (|2|) and (§) determine the behavior of the wave function in two distinct regions 
of space: the former in the interior of the cloud, the latter in the boundary region. For 
sufficiently large N these two regions are sufficiently extended to match each other. An 
example is shown in Fig. 1 for N = 10 5 . 

A third interesting region is the one at large distances beyond the boundary R where 
the system is very dilute and one can ignore the interaction term in Eq. (|I|). In this region 
the wave function can be written in the following way || 




(10) 



by 




(11) 
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The effects of the interatomic interactions enter here only through the value of the chemical 
potential. It is worth noticing that the case V ex t = would correspond to the asymptotic 
behavior of the order parameter for saturating systems in the absence of confining forces 



as happens, for example, outside the free surface of superfluid helium [10]. The coefficient 
of proportionality in (|1^) has been fixed in order to obtain the proper matching with the 
solution of Eqs.(|8|-[5p in the region of r such that R^> r — R^> d (see Eq.(PH|)). 

Let us apply the formalism discussed above to the simplest case of an isotropic harmonic 
trap: 

V ext (r) = \rnw 2 Ho r 2 . (13) 
For r < R, the Thomas-Fermi wave function (H) takes the form 



ipTF(r) 



(i? 2 -r^ 1/2 



(14) 



where we have used the expression \x = {l/2)muo 2 Ho R 2 for the chemical potential and intro- 
duced the harmonic oscillator length a HO = {Ti/muj HO ) l l 2 . The radius R is fixed by imposing 
the normalization of the wave function ( |14|) to the total number of particles: 

R 5 

N = ——— (15) 
15a< Q 

and increases very slowly with N. 

Near the boundary the wave function is instead given by Eq. (j|) where the thickness d, 
from Eq. ([/]), is 

,4 \ V3 

HO 

2R 



d =\^Tr) ■ (16) 



A similar result for the boundary thickness has been recently found by Baym and Pethick 
(see note 14 in Ref. ||). Taking large and negative values of £ as in Eq. flO]), means 
moving from the boundary to the interior of the cloud until (R — r) ^> d. In this region the 
asymptotic behavior (fllD holds and one obtains 



R(R - r) 
Arcaa 4 

HO 



1/2 



(17) 



This exactly coincides with ipTF given in Eq. fll4|) provided (R — r) <C R. In conclusion the 
wave function in the boundary region properly matches the Thomas-Fermi wave function 



(nj) for values of r satisfying the conditions 

d<&(R-r)<R. (18) 

For distances from the boundary less than d the Thomas-Fermi approximation (|T4] ) fails; 
vice versa, for distances comparable to the radius R, Eq. (||) becomes inadequate. 

Let us apply the above results to the calculation of the kinetic energy of the system. The 
integral (|3]) can be naturally divided into two parts: 

E kin = ^U W{r)\ 2 r 2 dr + J r ^ W(r)\ 2 r 2 dr\ (19) 

where the distance e > from the boundary R is chosen in such a way that the conditions 
(|18|), with (R — r) = e, are satisfied. This permits to evaluate the first term using the 



Thomas- Fermi approximation (p!^), and the second one using the solution ([8]) holding near 
the boundary. Clearly the sum of the two terms should not depend on the explicit value of 
e. 



(20) 



The first integral of Eq. (|T9|) is easily evaluated and becomes 



R-e c>3 







\ip'{r)\ 2 r 2 dr 
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where we have neglected corrections vanishing as e/R. 

For the second contribution, arising from the boundary region, we instead find the result 

/■+oo f?3 /"+oo 

/ W{r)\ 2 r 2 dr = —^ (<P'fdi (21) 
J R-e 4na^ HO a J-e/d 

where 0' = d<^jdl\. If the ratio e/d = (R — r)/d is sufficiently large (see condition ([18])) the 
integral in the right hand side is easily calculated and takes the value 

{<P') 2 di=-Aog^ + C (22) 



with 

C = - log (y/T+? + C) ^ ((0') 2 \/l+^) « = °- 176 (23) 
In Eq. (|23"D we have ignored corrections vanishing as d/e. Collecting the above results and 



using the explicit expression (|T^) for the boundary thickness d in terms of the oscillator 
length a HO , one finally finds the following result for the kinetic energy per particle: 

Ekin 5 Tl 



N 2mR 2 



log — +C 
V a 



"HO . 



2^ lo Hi^:J (24) 



where C = (7/4) log 2 — 2 + 3C and we have used expression ( |15|) for iV. Equation ( |24|) 
provides the correct asymptotic behavior of the kinetic energy in the limit of large iV where 
R 3> cl ho . This is confirmed by the comparison with the exact value of the kinetic energy 
obtained by solving numerically the Gross-Pitaevskii equation ([!]), as shown in Fig. 3. 

We conclude this paper by discussing an interesting application of the formalism to 
a Josephson-type effect. The physical idea is to consider a confining potential with two 
wells separated by a barrier. When the chemical potential in the two traps is different an 
oscillating flux of atoms is generated. Let us consider the simplest one-dimensional problem 
(extension to 3D will be the object of a future work) and let the external field V ext consist 
of two symmetric traps, trap 1 and trap 2, as shown schematically in Fig. 4. A difference 
between the chemical potentials Hi and \ii of the atoms in the two traps can be achieved, 
for example, by filling them with a different number of atoms. In order to obtain a first 
analytic result for the flux of atoms generated by the difference in the chemical potentials we 
will assume that the barrier between the two wells is high enough. In this case the overlap 
between the wave functions relative to the two traps, occurs only in the classically forbidden 
region where interaction effects can be ignored and one can safely use approximation ( |I~2"D for 
the wave function. Furthermore we will ignore the variation of fii and fi2 generated by the 
corresponding flux of particles. In ID the factor 1/r in the wave function (|1^) is absent and 
it is convenient to take the origin of axes at the symmetry point of the external potential 
(see Fig. 4). 



The Gross-Pitaevskii equation has two natural solutions in this case. The first one with 
chemical potential /ii, is localized in the trap 1. Its behavior in the classically forbidden 
region x > —L\ is given by 



1/2 



«*> = (iSfe) [(MVU*) -*)]■'' {-fwl'j^ - ^ /2dx 'j • < 25) 

where Xi is the distance between the center of trap 1 and the classical turning point, d\ 
is its boundary thickness [see Eq.(0)] and L\ is the distance between the classical turning 
point and the symmetry point of the external potential (see Fig. 4). 

The second solution with chemical potential \xi, is instead localized in the trap 2 and its 
behavior in the region x < L 2 is given by 



hX2 \ 1/2 A 



= \ydfa) [(2 m( K„W-^)]^ exp 1 1 lv '" (x,) - ^' 2dx ' 1 • (26) 

It is immediate to verify that the linear combination 

ip(x,t) =tp 1 (x)exp (-i^j-) +ij 2 (x)exp (-i^j (27) 

is solution of the time dependent Schrodinger equation. In fact the wave functions ipi an d 
■02 significantly overlap only in the classically forbidden region where non linear effects due 
to the interatomic potential are negligible. The current density 

/= ^ (W)|^*M) -nx,t)^(x,tfj (28) 

associated with the wave function ( p7|) can be easily calculated and takes the typical Joseph- 
son form 

I = Iq sm (29) 



with I = (h/rn^tyiip^—faipi)- Using the explicit results (p5|) and fl26|) for the wave functions 
ipi and ip2 and taking in ~ {1% — A* an d L\ ~ L2 = L in the evaluation of Iq, we find that 
the current Jo is uniform in the interval (-L, +L). Its explicit value is given by the useful 
result 
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hA 2 X 2 ( [2m f+ L r l1/9 , A 

As a consequence of the Josephson current the number of atoms in the two traps will 



oscillate in time according to the law [ITT 



d AT d AT T ■ {Pl-fr)* /Qi\ 

thereby providing the anticipated result for the flux of particles through the barrier sepa- 
rating the two traps. 

In conclusion we have obtained an explicit solution of the Gross-Pitaevskii equation 
near the classical turning point where the Thomas-Fermi approximation turns out to be 
completely inadequate. Using this solution we have been able to derive an analytic expression 
for the kinetic energy of the system holding for large values of N. We have finally discussed 
possible Josephson-type oscillations of atoms through the barrier separating two traps. 
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FIGURES 

FIG. 1. Condensate wave function for 10 5 atoms of 87 Rb (scattering length a = 5.29 x 1CP 7 
cm) in a spherical harmonic trap of length a HO = 1.22 x 10 -4 cm. Solid line: numerical solution 
of the Gross-Pitaevskii equation ([!]). Dot-dashed line: Thomas- Fermi approximation (^) (indis- 
tinguishable from the solid line in the inner part). Dashed line: surface profile obtained from the 
universal equation (||). 

FIG. 2. Solution of the universal equation (|9|). The two asymptotic limits (|Il| ) (dot-dashed 
line) and (|l0| ) (dashed line) are also shown. 

FIG. 3. Kinetic energy per particle, in units HcohOj f° r 87 Rb in a spherical harmonic trap as 
a function of the number of condensed atoms. Solid line: from the solution of the Gross-Pitaevskii 
equation (|l|). Dashed line: approximation (p4J). 

FIG. 4. Geometry of the double trap for the Josephson effect (see text). 
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